Collaborators: Shane Blowes, Jon Chase, Helmut Hillebrand, Michael Burrows, Amanda Bates, Uli Brose, Benoit Gauzens, Laura Antao Assistance: Katherine Lew, Josef Hauser
library(data.table) # for handling large datasets
library(ggplot2) # for some plotting
library(nlme) # for ME models
library(beanplot) # for beanplots
library(maps) # for map
library(gridExtra) # to combine ggplots together
library(grid) # to combine ggplots together
options(width=500) # turn off most text wrapping
# tell RStudio to use project root directory as the root for this notebook. Needed since we are storing code in a separate directory.
knitr::opts_knit$set(root.dir = rprojroot::find_rstudio_root_file())
# Turnover and covariates assembled by turnover_vs_temperature_prep.Rmd
trends <- fread('output/turnover_w_covariates.csv.gz')
# set realm order
trends[, REALM := factor(REALM, levels = c('Freshwater', 'Marine', 'Terrestrial'), ordered = FALSE)]
# realm that combined Terrestrial and Freshwater, for interacting with human impact
trends[, REALM2 := REALM]
levels(trends$REALM2) = list(TerrFresh = "Freshwater", TerrFresh = "Terrestrial", Marine = "Marine")
# group Marine invertebrates/plants in with All
trends[, taxa_mod2 := taxa_mod]
trends[taxa_mod == 'Marine invertebrates/plants', taxa_mod2 := 'All']
Log-transform some variables, then center and scale.
trends[, tempave.sc := scale(tempave)]
trends[, tempave_metab.sc := scale(tempave_metab)]
trends[, seas.sc := scale(seas)]
trends[, microclim.sc := scale(log(microclim))]
trends[, temptrend.sc := scale(temptrend)]
trends[, temptrend_abs.sc := scale(abs(temptrend))]
trends[, npp.sc := scale(log(npp))]
trends[, mass.sc := scale(log(mass_mean_weight))]
trends[, speed.sc := scale(log(speed_mean_weight+1))]
trends[, lifespan.sc := scale(log(lifespan_mean_weight))]
trends[, thermal_bias.sc := scale(thermal_bias)]
trends[, consumerfrac.sc := scale(consfrac)]
trends[, endothermfrac.sc := scale(endofrac)]
trends[, nspp.sc := scale(log(Nspp))]
trends[, human.sc := scale(human)]
# histograms to examine
cexmain = 0.6
par(mfrow = c(3,5))
invisible(trends[, hist(tempave.sc, main = 'Environmental temperature (°C)', cex.main = cexmain)])
invisible(trends[, hist(tempave_metab.sc, main = 'Metabolic temperature (°C)', cex.main = cexmain)])
invisible(trends[, hist(seas.sc, main = 'Seasonality (°C)', cex.main = cexmain)])
invisible(trends[, hist(microclim.sc, main = 'log Microclimates (°C)', cex.main = cexmain)])
invisible(trends[, hist(temptrend.sc, main = 'Temperature trend (°C/yr)', cex.main = cexmain)])
invisible(trends[, hist(temptrend_abs.sc, main = 'log abs(Temperature trend) (°C/yr)', cex.main = cexmain)])
invisible(trends[, hist(mass.sc, main = 'log Mass (g)', cex.main = cexmain)])
invisible(trends[, hist(speed.sc, main = 'log Speed (km/hr)', cex.main = cexmain)])
invisible(trends[, hist(lifespan.sc, main = 'log Lifespan (yr)', cex.main = cexmain)])
invisible(trends[, hist(consumerfrac.sc, main = 'Consumers (fraction)', cex.main = cexmain)])
invisible(trends[, hist(endothermfrac.sc, main = 'Endotherms (fraction)', cex.main = cexmain)])
invisible(trends[, hist(nspp.sc, main = 'log Species richness', cex.main = cexmain)])
invisible(trends[, hist(thermal_bias.sc, main = 'Thermal bias (°C)', cex.main = cexmain)])
invisible(trends[, hist(npp.sc, main = 'log Net primary productivity', cex.main = cexmain)])
invisible(trends[, hist(human.sc, main = 'Human impact score', cex.main = cexmain)])
panel.cor <- function(x, y, digits = 2, prefix = "", cex.cor, ...)
{
usr <- par("usr"); on.exit(par(usr))
par(usr = c(0, 1, 0, 1))
r <- cor(x, y, use = 'pairwise.complete.obs')
txt <- format(c(r, 0.123456789), digits = digits)[1]
txt <- paste0(prefix, txt)
if(missing(cex.cor)) cex.cor <- 0.8/strwidth(txt)
text(0.5, 0.5, txt) #, cex = cex.cor * r)
}
pairs(formula = ~ REALM + tempave.sc + tempave_metab.sc + seas.sc + microclim.sc + temptrend.sc + temptrend_abs.sc + mass.sc + speed.sc + lifespan.sc + consumerfrac.sc + endothermfrac.sc + nspp.sc + thermal_bias.sc + npp.sc + human.sc, data = trends, gap = 1/10, cex = 0.2, col = '#00000022', lower.panel = panel.cor)
Mass and lifespan look tightly correlated, but r only 0.56…? Tempave_metab and lifespan don’t look tightly correlated, but r= -0.81 Tempave_metab and speed don’t look tightly correlated, but r= -0.83 Lifespan and speed don’t look tightly correlated, but r = 0.73
Just turnover
cat('Overall # time-series: ', nrow(trends), '\n')
Overall # time-series: 53013
cat('# studies: ', trends[, length(unique(STUDY_ID))], '\n')
# studies: 332
cat('Data points: ', trends[, sum(nyrBT)], '\n')
Data points: 293973
trends[, table(REALM)]
REALM
Freshwater Marine Terrestrial
1025 48647 3341
trends[, table(taxa_mod)]
taxa_mod
All Amphibians Benthos Birds
1705 379 4679 13741
Fish Invertebrates Mammals Marine invertebrates/plants
28473 2996 525 206
Plant Reptiles
305 4
trends[, table(taxa_mod, REALM)]
REALM
taxa_mod Freshwater Marine Terrestrial
All 0 1702 3
Amphibians 2 0 377
Benthos 0 4679 0
Birds 0 11099 2642
Fish 1006 27467 0
Invertebrates 15 2901 80
Mammals 0 478 47
Marine invertebrates/plants 0 206 0
Plant 1 115 189
Reptiles 1 0 3
With all covariates
# the cases we can compare
apply(trends[, .(Jtutrend, REALM, tempave.sc, tempave_metab.sc, seas.sc, microclim.sc, temptrend.sc, mass.sc, speed.sc, lifespan.sc, consumerfrac.sc, endothermfrac.sc, nspp.sc, thermal_bias.sc, npp.sc, human.sc)], MARGIN = 2, FUN = function(x) sum(!is.na(x)))
Jtutrend REALM tempave.sc tempave_metab.sc seas.sc microclim.sc
53013 53013 49916 49916 49916 51834
temptrend.sc mass.sc speed.sc lifespan.sc consumerfrac.sc endothermfrac.sc
49916 52820 52689 51540 47534 53013
nspp.sc thermal_bias.sc npp.sc human.sc
53013 49371 52863 53013
i <- trends[, complete.cases(Jtutrend, temptrend.sc, tempave_metab.sc, REALM, seas.sc, microclim.sc, npp.sc, mass.sc, speed.sc, lifespan.sc, consumerfrac.sc, thermal_bias.sc)]
cat('Overall # time-series: ', sum(i), '\n')
Overall # time-series: 43585
cat('# studies: ', trends[i, length(unique(STUDY_ID))], '\n')
# studies: 250
cat('Data points: ', trends[i, sum(nyrBT)], '\n')
Data points: 222824
trends[i, table(REALM)]
REALM
Freshwater Marine Terrestrial
1008 39735 2842
trends[i, table(taxa_mod)]
taxa_mod
All Amphibians Benthos Birds Fish Invertebrates Mammals Plant
521 12 590 11803 27372 2567 518 200
Reptiles
2
trends[i, table(taxa_mod, REALM)]
REALM
taxa_mod Freshwater Marine Terrestrial
All 0 520 1
Amphibians 2 0 10
Benthos 0 590 0
Birds 0 9221 2582
Fish 993 26379 0
Invertebrates 12 2484 71
Mammals 0 477 41
Plant 1 64 135
Reptiles 0 0 2
Try combinations of
And choose the one with lowest AIC (not run: takes a long time)
# fit models for variance structure
fixed <- formula(Jtutrend ~ REALM + tempave_metab.sc + seas.sc + microclim.sc + npp.sc + temptrend_abs.sc +
mass.sc + speed.sc + lifespan.sc + consumerfrac.sc + thermal_bias.sc)
i <- trends[, complete.cases(Jtutrend, REALM, tempave_metab.sc, seas.sc, microclim.sc, npp.sc, temptrend_abs.sc,
mass.sc, speed.sc, lifespan.sc, consumerfrac.sc, thermal_bias.sc)]
mods <- vector('list', 0)
mods[[1]] <- gls(fixed, data = trends[i,])
mods[[2]] <- gls(fixed, data = trends[i,], weights = varPower(-0.5, ~nyrBT))
mods[[3]] <- gls(fixed, data = trends[i,], weights = varPower(0.5, ~ abs(temptrend)))
mods[[4]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2, control = lmeControl(opt = "optim"))
mods[[5]] <- lme(fixed, data = trends[i,], random = ~1|STUDY_ID, control = lmeControl(opt = "optim"))
mods[[6]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2/STUDY_ID, control = lmeControl(opt = "optim"))
mods[[7]] <- lme(fixed, data = trends[i,], random = ~1|STUDY_ID/rarefyID, control = lmeControl(opt = "optim"))
mods[[8]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2/STUDY_ID/rarefyID, control = lmeControl(opt = "optim"))
mods[[9]] <- lme(fixed, data = trends[i,], random = ~temptrend_abs.sc | taxa_mod)
mods[[10]] <- lme(fixed, data = trends[i,], random = ~temptrend_abs.sc | STUDY_ID)
mods[[11]] <- lme(fixed, data = trends[i,], random = ~temptrend_abs.sc | taxa_mod2/STUDY_ID, control = lmeControl(opt = "optim"))
mods[[12]] <- lme(fixed, data = trends[i,], random = list(STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1)) # includes overdispersion. new formula so that random slope is only for study level (not enough data to extend to rarefyID).
mods[[13]] <- lme(fixed, data = trends[i,], random = list(taxa_mod2 = ~ temptrend_abs.sc, STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1)) # 30+ min to fit
mods[[14]] <- lme(fixed, data = trends[i,], random = ~1|STUDY_ID, weights = varPower(-0.5, ~nyrBT))
mods[[15]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2, weights = varPower(-0.5, ~nyrBT))
mods[[16]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2/STUDY_ID, weights = varPower(-0.5, ~nyrBT))
mods[[17]] <- lme(fixed, data = trends[i,], random = ~1|STUDY_ID/rarefyID, weights = varPower(-0.5, ~nyrBT))
mods[[18]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2/STUDY_ID/rarefyID, weights = varPower(-0.5, ~nyrBT))
mods[[19]] <- lme(fixed, data = trends[i,], random = ~temptrend_abs.sc|STUDY_ID, weights = varPower(-0.5, ~nyrBT))
mods[[20]] <- lme(fixed, data = trends[i,], random = list(STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1), weights = varPower(-0.5, ~nyrBT))
mods[[21]] <- lme(fixed, data = trends[i,], random = list(taxa_mod2 = ~ temptrend_abs.sc, STUDY_ID = ~ 1), weights = varPower(-0.5, ~nyrBT))
mods[[22]] <- lme(fixed, data = trends[i,], random = list(taxa_mod2 = ~ temptrend_abs.sc, STUDY_ID = ~ 1, rarefyID = ~1), weights = varPower(-0.5, ~nyrBT))
mods[[23]] <- lme(fixed, data = trends[i,], random = list(taxa_mod2 = ~ temptrend_abs.sc, STUDY_ID = ~ temptrend_abs.sc), weights = varPower(-0.5, ~nyrBT)) # singular precision warning with lmeControl(opt = 'optim') and convergence error without
mods[[24]] <- lme(fixed, data = trends[i,], random = list(taxa_mod2 = ~ temptrend_abs.sc, STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1), weights = varPower(-0.5, ~nyrBT)) # singular precision warning with lmeControl(opt = 'optim') and convergence error without
mods[[25]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2, weights = varPower(-0.5, ~abs(temptrend)))
mods[[26]] <- lme(fixed, data = trends[i,], random = ~1|STUDY_ID, weights = varPower(-0.5, ~abs(temptrend)))
mods[[27]] <- lme(fixed, data = trends[i,], random = ~1|STUDY_ID/rarefyID, weights = varPower(-0.5, ~abs(temptrend)))
mods[[28]] <- lme(fixed, data = trends[i,], random = ~1|taxa_mod2/STUDY_ID/rarefyID, weights = varPower(-0.5, ~abs(temptrend)))
mods[[29]] <- lme(fixed, data = trends[i,], random = ~temptrend_abs.sc|STUDY_ID, weights = varPower(-0.5, ~abs(temptrend)))
mods[[30]] <- lme(fixed, data = trends[i,], random = ~temptrend_abs.sc|taxa_mod2/STUDY_ID, weights = varPower(-0.5, ~abs(temptrend)), control = lmeControl(opt = "optim"))
mods[[31]] <- lme(fixed, data = trends[i,], random = list(STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1), weights = varPower(-0.5, ~abs(temptrend)))
mods[[32]] <- lme(fixed, data = trends[i,], random = list(taxa_mod2 = ~ temptrend_abs.sc, STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1), weights = varPower(-0.5, ~abs(temptrend)), control = lmeControl(opt = "optim")) # singular precision warning
aics <- sapply(mods, AIC)
minaics <- aics - min(aics)
minaics
which.min(aics)
Chooses the random slopes (temptrend_abs) & intercepts for STUDY_ID, overdispersion, and variance scaled to number of years. We haven’t dealt with potential testing on the boundary issues here yet.
world <- map_data('world')
ggplot(world, aes(x = long, y = lat, group = group)) +
geom_polygon(fill = 'lightgray', color = 'white') +
geom_point(data = trends, aes(rarefyID_x, rarefyID_y, group = REALM, color = REALM), size = 0.5, alpha = 0.4) +
scale_color_brewer(palette="Set1", name = 'Realm') +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.key=element_blank(),
axis.text=element_text(size=16),
axis.title=element_text(size=20)) +
labs(x = 'Longitude (°)', y = 'Latitude (°)')
Mostly northern hemisphere, but spread all over. No so much in Africa or much of Asia.
Lines are ggplot smoother fits by realm.
Strong trends with temperature change, but trends are pretty symmetric around no trend in temperature, which implies warming or cooling drives similar degree of community turnover. Some indication of less turnover for larger organisms (mass) Higher turnover on land with higher seasonality? More turnover for shorter-lived organisms? No really clear differences among realms.
Average rates of turnover
trends[abs(temptrend) >= 0.5, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # turnover per year for locations changing temperature
trends[abs(temptrend) < 0.1, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # not changing temperature
trends[temptrend >= 0.5, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # warming
trends[temptrend <= -0.5, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # cooling
trends[abs(temptrend) >= 0.5 & abs(rarefyID_y) < 35, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # tropics and sub-tropics
trends[abs(temptrend) >= 0.5 & abs(rarefyID_y) >= 35 & abs(rarefyID_y) < 66.56339, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # temperate
trends[abs(temptrend) >= 0.5 & abs(rarefyID_y) >= 66.56339, .(mean(Jtutrend), sd(Jtutrend)/sqrt(.N))] # arctic
i <- trends[, !duplicated(rarefyID)]; sum(i)
[1] 53013
par(mfrow=c(5,3))
beanplot(rarefyID_y ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Latitude (degN)', ll = 0.05)
beanplot(tempave ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Temperature (degC)', ll = 0.05)
beanplot(tempave_metab ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Metabolic Temperature (degC)', ll = 0.05, bw = 'nrd0') # nrd0 bandwidth to calculation gap
beanplot(seas ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Seasonality (degC)', ll = 0.05)
beanplot(microclim ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Microclimates (degC)', ll = 0.05)
log="y" selected
beanplot(temptrend ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Temperature trend (degC/yr)', ll = 0.05)
beanplot(mass_mean_weight ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Mass (g)', ll = 0.05, log = 'y')
beanplot(speed_mean_weight +1 ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Speed (km/hr)', ll = 0.05, log = 'y')
beanplot(lifespan_mean_weight ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Lifespan (yr)', ll = 0.05, log = 'y')
#beanplot(consfrac ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Consumers (fraction)', ll = 0.05, log = '') # too sparse
#beanplot(endofrac ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Endotherms (fraction)', ll = 0.05, log = '') # too sparse
beanplot(Nspp ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Number of species', ll = 0.05, log = 'y')
beanplot(thermal_bias ~ REALM, data = trends[i & !is.na(thermal_bias),], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Thermal bias (degC)', ll = 0.05)
beanplot(npp ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'NPP', ll = 0.05)
log="y" selected
beanplot(human ~ REALM, data = trends[i,], what = c(1,1,1,1), col = c("#CAB2D6", "#33A02C", "#B2DF8A"), border = "#CAB2D6", ylab = 'Human impact score', ll = 0.05)
Marine are in generally warmer locations (seawater doesn’t freeze) Marine have much lower seasonality. Marine and freshwater have some very small masses (plankton), but much of dataset is similar to terrestrial. Marine has a lot of slow, crawling organisms, but land has plants. Land also has birds (fast).
i <- trends[, complete.cases(Jtutrend, REALM, temptrend)]
randef <- list(STUDY_ID = ~ abs(temptrend), rarefyID = ~1)
varef <- varPower(-0.5, ~nyrBT)
if(file.exists('temp/modonlyTtrend.rds')){
modonlyTtrend <- readRDS('temp/modonlyTtrend.rds')
} else {
modonlyTtrend <- lme(Jtutrend ~ abs(temptrend)*REALM,
random = randef, weights = varef, data = trends[i,], method = 'REML')
saveRDS(modonlyTtrend, file = 'temp/modonlyTtrend.rds')
}
i2 <- trends[, complete.cases(Jbetatrend, REALM, temptrend)]
if(file.exists('temp/modonlyTtrendJbeta.rds')){
modonlyTtrendJbeta <- readRDS('temp/modonlyTtrendJbeta.rds')
} else {
modonlyTtrendJbeta <- lme(Jbetatrend ~ abs(temptrend)*REALM,
random = randef, weights = varef, data = trends[i2,], method = 'REML',
control=lmeControl(msMaxIter = 100, maxIter = 100))
saveRDS(modonlyTtrendJbeta, file = 'temp/modonlyTtrendJbeta.rds')
}
i3 <- trends[, complete.cases(Horntrend, REALM, temptrend)]
if(file.exists('temp/modonlyTtrendHorn.rds')){
modonlyTtrendHorn <- readRDS('temp/modonlyTtrendHorn.rds')
} else {
modonlyTtrendHorn <- lme(Horntrend ~ abs(temptrend)*REALM,
random = randef, weights = varef, data = trends[i3,], method = 'REML')
saveRDS(modonlyTtrendHorn, file = 'temp/modonlyTtrendHorn.rds')
}
summary(modonlyTtrend)
Linear mixed-effects model fit by REML
Data: trends[i, ]
Random effects:
Formula: ~abs(temptrend) | STUDY_ID
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.04375688 (Intr)
abs(temptrend) 0.29351811 -0.123
Formula: ~1 | rarefyID %in% STUDY_ID
(Intercept) Residual
StdDev: 0.001792571 0.2893301
Variance function:
Structure: Power of variance covariate
Formula: ~nyrBT
Parameter estimates:
power
-1.223838
Fixed effects: Jtutrend ~ abs(temptrend) * REALM
Correlation:
(Intr) abs(t) REALMM REALMT a():REALMM
abs(temptrend) -0.424
REALMMarine -0.934 0.396
REALMTerrestrial -0.923 0.391 0.862
abs(temptrend):REALMMarine 0.407 -0.960 -0.407 -0.375
abs(temptrend):REALMTerrestrial 0.387 -0.913 -0.362 -0.433 0.876
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-8.38476416 -0.29781039 0.08738927 0.54929160 6.52538535
Number of Observations: 49916
Number of Groups:
STUDY_ID rarefyID %in% STUDY_ID
317 49916
summary(modonlyTtrendJbeta)
Linear mixed-effects model fit by REML
Data: trends[i2, ]
Random effects:
Formula: ~abs(temptrend) | STUDY_ID
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.05691594 (Intr)
abs(temptrend) 0.34739161 -0.16
Formula: ~1 | rarefyID %in% STUDY_ID
(Intercept) Residual
StdDev: 1.28009e-06 0.3525846
Variance function:
Structure: Power of variance covariate
Formula: ~nyrBT
Parameter estimates:
power
-1.394679
Fixed effects: Jbetatrend ~ abs(temptrend) * REALM
Correlation:
(Intr) abs(t) REALMM REALMT a():REALMM
abs(temptrend) -0.428
REALMMarine -0.933 0.399
REALMTerrestrial -0.926 0.396 0.864
abs(temptrend):REALMMarine 0.411 -0.961 -0.411 -0.380
abs(temptrend):REALMTerrestrial 0.390 -0.913 -0.364 -0.437 0.877
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.3578768 -0.1587203 0.2033006 0.6705745 6.5270460
Number of Observations: 49916
Number of Groups:
STUDY_ID rarefyID %in% STUDY_ID
317 49916
summary(modonlyTtrendHorn)
Linear mixed-effects model fit by REML
Data: trends[i3, ]
Random effects:
Formula: ~abs(temptrend) | STUDY_ID
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.0549979 (Intr)
abs(temptrend) 0.3228314 -0.095
Formula: ~1 | rarefyID %in% STUDY_ID
(Intercept) Residual
StdDev: 0.007989212 0.3034899
Variance function:
Structure: Power of variance covariate
Formula: ~nyrBT
Parameter estimates:
power
-1.217439
Fixed effects: Horntrend ~ abs(temptrend) * REALM
Correlation:
(Intr) abs(t) REALMM REALMT a():REALMM
abs(temptrend) -0.399
REALMMarine -0.921 0.368
REALMTerrestrial -0.919 0.367 0.847
abs(temptrend):REALMMarine 0.381 -0.954 -0.378 -0.350
abs(temptrend):REALMTerrestrial 0.363 -0.908 -0.334 -0.408 0.867
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-6.21250285 -0.32193926 0.06521414 0.53688251 6.03752738
Number of Observations: 48800
Number of Groups:
STUDY_ID rarefyID %in% STUDY_ID
276 48800
# make table of coefficients
coefs <- as.data.frame(summary(modonlyTtrend)$tTable)
coefs2 <- as.data.frame(summary(modonlyTtrendJbeta)$tTable)
coefs3 <- as.data.frame(summary(modonlyTtrendHorn)$tTable)
coefs$mod <- 'Jtu'
coefs2$mod <- 'Jbeta'
coefs3$mod <- 'Horn'
rows1 <- which(grepl('temptrend', rownames(coefs))) # extract temperature effect
cols <- c('Value', 'Std.Error', 'mod')
allcoefs <- rbind(coefs[rows1, cols], coefs2[rows1, cols], coefs3[rows1, cols])
allcoefs$Value[grepl('REALMMarine', rownames(allcoefs))] <-
allcoefs$Value[grepl('REALMMarine', rownames(allcoefs))] +
allcoefs$Value[!grepl('REALM', rownames(allcoefs))] # add intercept to marine effects
allcoefs$Value[grepl('REALMTerrestrial', rownames(allcoefs))] <-
allcoefs$Value[grepl('REALMTerrestrial', rownames(allcoefs))] +
allcoefs$Value[!grepl('REALM', rownames(allcoefs))] # add intercept to terrestrial effects
allcoefs$lCI <- allcoefs$Value - allcoefs$Std.Error # lower confidence interval
allcoefs$uCI <- allcoefs$Value + allcoefs$Std.Error
allcoefs$y <- c(3, 2, 1, 2.9, 1.9, 0.9, 2.8, 1.8, 0.8) # y-values
allcoefs$col <- c(rep('black', 3), rep('light grey', 3), rep('dark grey', 3))
allcoefs$realm <- rep(c('Freshwater', 'Marine', 'Terrestrial'), 3)
par(las = 1, mai = c(0.8, 2, 0.1, 0.1))
plot(0,0, col = 'white', xlim=c(-0.02, 0.85), ylim = c(0.7,3),
yaxt='n', xlab = 'Turnover per |°C/yr|', ylab ='')
axis(2, at = 3:1, labels = c('Freshwater', 'Marine', 'Terrestrial'), cex.axis = 0.7)
abline(v = 0, col = 'grey')
for(i in 1:nrow(allcoefs)){
with(allcoefs[i, ], points(Value, y, pch = 16, col = col))
with(allcoefs[i, ], lines(x = c(lCI, uCI), y = c(y, y), col = col))
}
legend('bottomright', col = c('black', 'dark grey', 'light grey'), lwd = 1, pch = 16,
legend = c('Jaccard turnover', 'Jaccard total', 'Horn-Morisita'))
Scatterplot, violin plots, and coefficient plots all together
# on macbook: fig.width=3, fig.height=2.375, fig.retina=3, out.width=3, out.height=2.375
# on external monitor:
trends[temptrend <= -0.5, temptrendtext := 'Cooling']
trends[abs(temptrend) <= 0.1, temptrendtext := 'Stable']
trends[temptrend >= 0.5, temptrendtext := 'Warming']
trends[abs(rarefyID_y) < 35, latzone := 'Subtropics']
trends[abs(rarefyID_y) >= 35 & abs(rarefyID_x) < 66.56339, latzone := 'Temperate']
trends[abs(rarefyID_y) >= 66.56339, latzone := 'Polar']
trends[, latzone := factor(latzone, levels = c('Subtropics', 'Temperate', 'Polar'))]
p1 <- ggplot(trends, aes(temptrend, Jtutrend, color = REALM, fill = REALM, size = nyrBT)) +
geom_point(na.rm = TRUE, shape = 16, alpha = 0.1) +
geom_smooth(data=subset(trends, abs(temptrend) < 0.75), method = 'gam', formula = y ~ s(x, bs = "cs"),
na.rm = TRUE) +
scale_color_brewer(palette="Set1", name = 'Realm') +
scale_fill_brewer(palette="Set1", name = 'Realm') +
labs(x = 'Temperature trend (°C/yr)', y = 'Jaccard turnover', tag = 'A') +
scale_size_continuous(range = c(1, 8), breaks = c(2, 5, 20)) +
guides(size = guide_legend(title = 'Years',
override.aes = list(linetype=0, fill = NA, alpha = 1))) +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.key=element_blank(),
legend.key.size = unit(0.5,"line"),
axis.text=element_text(size=8),
axis.title=element_text(size=10))
p2 <- ggplot(trends[!is.na(temptrendtext), ], aes(temptrendtext, Jtutrend)) +
geom_violin(draw_quantiles = c(0.25, 0.5, 0.75), fill = 'grey') +
labs(x = '', y = 'Jaccard turnover', tag = 'B') +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.key=element_blank(),
axis.text=element_text(size=8),
axis.title=element_text(size=10))
p3 <- ggplot(trends[abs(temptrend) >= 0.5 & !is.na(latzone), ], aes(latzone, Jtutrend)) +
geom_violin(draw_quantiles = c(0.25, 0.5, 0.75), fill = 'grey') +
labs(x = '', y = '', tag = 'C') +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.key=element_blank(),
axis.text=element_text(size=7),
axis.title=element_text(size=10))
p4 <- ggplot(allcoefs, aes(Value, y, group = mod, color = mod)) +
geom_errorbarh(aes(xmin = lCI, xmax = uCI, height = 0)) +
geom_point() +
labs(x = expression(atop('Temperature change effect', '(Turnover '~degree*'C'^'-1'*')')), y = '', tag = 'D') +
scale_color_grey() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.position='none',
axis.text=element_text(size=7),
axis.title=element_text(size=7)) +
coord_cartesian(xlim =c(0, 1)) +
scale_y_continuous(name = '', breaks = c(1, 2, 3), labels = c('Terrestrial', 'Marine', 'Freshwater'))
grid.arrange(p1, p2, p3, p4, ncol = 3, layout_matrix = rbind(c(1,1,1), c(2,3,4)),
heights=c(unit(0.66, "npc"), unit(0.34, "npc")))
Try static covariates plus interactions of abs temperature trend with each covariate:
Except for thermal bias: interact with temperature trend (not abs)
summary(modTfull1)
Linear mixed-effects model fit by REML
Data: trends[i, ]
Random effects:
Formula: ~temptrend_abs.sc | STUDY_ID
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.04834247 (Intr)
temptrend_abs.sc 0.03064449 0.335
Formula: ~1 | rarefyID %in% STUDY_ID
(Intercept) Residual
StdDev: 0.0008908894 0.2861661
Variance function:
Structure: Power of variance covariate
Formula: ~nyrBT
Parameter estimates:
power
-1.235108
Fixed effects: Jtutrend ~ temptrend_abs.sc * REALM + temptrend_abs.sc * tempave.sc + temptrend_abs.sc * tempave_metab.sc + temptrend_abs.sc * seas.sc + temptrend_abs.sc * microclim.sc + temptrend_abs.sc * mass.sc + temptrend_abs.sc * speed.sc + temptrend_abs.sc * lifespan.sc + temptrend_abs.sc * consumerfrac.sc + temptrend_abs.sc * endothermfrac.sc + temptrend_abs.sc * nspp.sc + temptrend.sc * thermal_bias.sc + temptrend_abs.sc * npp.sc + temptrend_abs.sc * human.sc:REALM2
Correlation:
(Intr) tmpt_. REALMM REALMT tmpv.s tmpv_. ses.sc mcrcl. mss.sc spd.sc lfspn. cnsmr. endth. nspp.s tmptr.
temptrend_abs.sc 0.344
REALMMarine -0.917 -0.314
REALMTerrestrial -0.878 -0.314 0.802
tempave.sc 0.002 -0.008 -0.015 0.002
tempave_metab.sc 0.047 0.030 -0.035 -0.050 -0.338
seas.sc -0.036 -0.019 0.051 -0.020 0.181 0.007
microclim.sc -0.012 -0.019 0.015 -0.016 -0.044 0.133 0.135
mass.sc 0.012 -0.012 0.000 0.003 -0.042 -0.515 -0.052 -0.027
speed.sc 0.019 0.002 -0.032 -0.032 -0.011 0.121 -0.059 -0.015 -0.090
lifespan.sc 0.046 0.047 -0.028 -0.048 -0.023 0.713 0.075 0.031 -0.813 0.190
consumerfrac.sc -0.072 -0.057 0.081 0.342 -0.001 -0.109 0.002 0.004 0.060 -0.164 -0.131
endothermfrac.sc 0.172 0.080 -0.091 -0.347 0.134 -0.234 0.004 -0.044 0.011 0.085 -0.024 -0.415
nspp.sc -0.013 0.021 -0.010 -0.023 -0.013 0.032 -0.033 -0.024 -0.169 -0.062 0.138 0.012 0.053
temptrend.sc -0.008 -0.015 0.008 0.008 -0.087 0.004 -0.034 0.018 0.017 -0.004 -0.014 0.005 -0.007 -0.019
thermal_bias.sc 0.029 -0.011 -0.045 -0.008 0.671 -0.024 -0.113 -0.112 0.033 0.028 -0.059 -0.002 -0.012 -0.069 0.000
npp.sc -0.007 -0.033 -0.001 0.008 -0.221 0.176 -0.229 -0.291 0.004 -0.012 0.044 -0.043 -0.070 -0.119 0.022
temptrend_abs.sc:REALMMarine -0.327 -0.952 0.337 0.296 0.006 -0.021 0.031 0.020 0.014 -0.006 -0.032 0.060 -0.055 -0.026 0.013
temptrend_abs.sc:REALMTerrestrial -0.297 -0.852 0.269 0.369 0.014 -0.048 -0.051 -0.025 0.015 -0.010 -0.040 0.151 -0.156 -0.041 0.010
temptrend_abs.sc:tempave.sc 0.001 -0.008 -0.004 0.013 0.489 -0.305 0.006 -0.082 -0.083 -0.040 0.014 0.026 0.151 0.021 -0.085
temptrend_abs.sc:tempave_metab.sc 0.014 0.050 -0.008 -0.024 -0.183 0.668 0.078 0.126 -0.299 0.033 0.434 -0.077 -0.191 0.030 0.013
temptrend_abs.sc:seas.sc -0.013 -0.038 0.018 -0.027 0.109 0.079 0.648 0.080 -0.045 -0.018 0.036 0.006 -0.034 -0.015 -0.005
temptrend_abs.sc:microclim.sc -0.012 -0.030 0.009 -0.012 -0.042 0.131 0.086 0.768 -0.033 0.010 0.013 -0.001 -0.056 0.013 0.029
temptrend_abs.sc:mass.sc -0.006 -0.014 0.008 0.009 -0.085 -0.330 -0.051 -0.028 0.605 0.007 -0.528 0.031 0.001 -0.105 0.020
temptrend_abs.sc:speed.sc 0.001 0.010 -0.003 -0.010 -0.025 0.015 -0.007 0.007 0.018 0.286 0.005 -0.057 0.039 -0.040 -0.010
temptrend_abs.sc:lifespan.sc 0.026 0.071 -0.013 -0.023 0.026 0.465 0.045 0.014 -0.515 0.021 0.645 -0.076 -0.034 0.112 -0.010
temptrend_abs.sc:consumerfrac.sc -0.054 -0.094 0.054 0.135 0.032 -0.122 -0.011 -0.002 0.039 -0.061 -0.111 0.360 -0.118 -0.002 0.003
temptrend_abs.sc:endothermfrac.sc 0.072 0.164 -0.046 -0.145 0.140 -0.270 -0.050 -0.077 -0.002 0.023 -0.035 -0.135 0.443 0.054 -0.018
temptrend_abs.sc:nspp.sc 0.012 0.036 -0.019 -0.027 -0.005 0.066 -0.030 0.010 -0.110 -0.029 0.125 -0.001 0.021 0.583 -0.029
temptrend.sc:thermal_bias.sc 0.006 0.006 -0.007 -0.004 0.076 0.038 0.031 -0.054 -0.010 0.013 0.006 -0.004 -0.013 -0.023 -0.354
thrm_. npp.sc t_.:REALMM t_.:REALMT tmptrnd_bs.sc:t. t_.:_. tmptrnd_bs.sc:ss. tmptrnd_bs.sc:mc.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc -0.085
temptrend_abs.sc:REALMMarine 0.009 0.019
temptrend_abs.sc:REALMTerrestrial 0.003 0.023 0.801
temptrend_abs.sc:tempave.sc 0.049 -0.066 0.006 0.027
temptrend_abs.sc:tempave_metab.sc 0.056 0.072 -0.038 -0.067 -0.502
temptrend_abs.sc:seas.sc 0.079 -0.136 0.054 -0.064 0.137 0.073
temptrend_abs.sc:microclim.sc -0.003 -0.256 0.034 -0.026 -0.060 0.106 0.113
temptrend_abs.sc:mass.sc -0.011 0.019 0.020 0.023 -0.117 -0.479 -0.047 -0.018
temptrend_abs.sc:speed.sc 0.002 0.004 -0.026 -0.031 -0.043 0.063 -0.016 -0.003
temptrend_abs.sc:lifespan.sc 0.004 0.020 -0.052 -0.064 0.028 0.667 0.042 0.010
temptrend_abs.sc:consumerfrac.sc -0.003 -0.034 0.103 0.309 0.046 -0.164 -0.026 0.000
temptrend_abs.sc:endothermfrac.sc -0.035 -0.033 -0.111 -0.365 0.369 -0.413 -0.048 -0.065
temptrend_abs.sc:nspp.sc -0.032 -0.037 -0.051 -0.071 0.034 0.040 -0.050 -0.022
temptrend.sc:thermal_bias.sc -0.023 -0.008 -0.007 -0.004 0.044 0.046 0.047 -0.040
tmptrnd_bs.sc:ms. tmptrnd_bs.sc:sp. tmptrnd_bs.sc:l. tmptrnd_bs.sc:c. tmptrnd_bs.sc:nd. tmptrnd_bs.sc:ns.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc -0.026
temptrend_abs.sc:lifespan.sc -0.831 0.093
temptrend_abs.sc:consumerfrac.sc 0.061 -0.188 -0.175
temptrend_abs.sc:endothermfrac.sc -0.018 0.126 -0.032 -0.292
temptrend_abs.sc:nspp.sc -0.208 -0.083 0.201 0.001 0.103
temptrend.sc:thermal_bias.sc -0.016 0.028 0.002 -0.003 -0.031 -0.029
tm.:_. tmptrnd_bs.sc:np. h.:REALM2T h.:REALM2M t_.:.:REALM2T
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc
temptrend_abs.sc:consumerfrac.sc
temptrend_abs.sc:endothermfrac.sc
temptrend_abs.sc:nspp.sc
temptrend.sc:thermal_bias.sc
[ reached getOption("max.print") -- omitted 5 rows ]
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-8.34115139 -0.34664125 0.05370777 0.53577626 6.61354550
Number of Observations: 43585
Number of Groups:
STUDY_ID rarefyID %in% STUDY_ID
250 43585
coefs <- summary(modTfull1)$tTable
par(las = 1, mai = c(0.5, 3, 0.1, 0.1))
rows1 <- which(!grepl('Intercept', rownames(coefs)))
plot(0,0, col = 'white', xlim=c(-0.05, 0.08), ylim = c(1,length(rows1)), yaxt='n', xlab = '', ylab ='')
axis(2, at = length(rows1):1, labels = rownames(coefs)[rows1], cex.axis = 0.7)
abline(v = 0, col = 'grey')
for(i in 1:length(rows1)){
x = coefs[rows1[i], 1]
se = coefs[rows1[i], 2]
points(x, length(rows1) + 1 - i, pch = 16)
lines(x = c(x-se, x+se), y = c(length(rows1) + 1 - i, length(rows1) + 1 - i))
}
resids <- resid(modTfull1)
preds <- getData(modTfull1)
col = '#00000033'
cex = 0.5
par(mfrow = c(4,4))
boxplot(resids ~ preds$REALM, cex = cex, col = col)
plot(preds$temptrend_abs.sc, resids, cex = cex, col = col)
plot(preds$temptrend.sc, resids, cex = cex, col = col)
plot(preds$tempave.sc, resids, cex = cex, col = col)
plot(preds$tempave_metab.sc, resids, cex = cex, col = col)
plot(preds$seas.sc, resids, cex = cex, col = col)
plot(preds$microclim.sc, resids, cex = cex, col = col)
plot(preds$mass.sc, resids, cex = cex, col = col)
plot(preds$speed.sc, resids, cex = cex, col = col)
plot(preds$lifespan.sc, resids, cex = cex, col = col)
plot(preds$consumerfrac.sc, resids, cex = cex, col = col)
plot(preds$endothermfrac.sc, resids, cex = cex, col = col)
plot(preds$nspp.sc, resids, cex = cex, col = col)
plot(preds$thermal_bias.sc, resids, cex = cex, col = col)
plot(preds$npp.sc, resids, cex = cex, col = col)
plot(preds$human.sc, resids, cex = cex, col = col)
i2 <- trends[, complete.cases(Jbetatrend, REALM, tempave.sc, tempave_metab.sc, seas.sc, microclim.sc,
temptrend.sc, temptrend_abs.sc, mass.sc, speed.sc, lifespan.sc,
consumerfrac.sc, endothermfrac.sc, nspp.sc, thermal_bias.sc, npp.sc, human.sc)]
i3 <- trends[, complete.cases(Horntrend, REALM, tempave.sc, tempave_metab.sc, seas.sc, microclim.sc,
temptrend.sc, temptrend_abs.sc, mass.sc, speed.sc, lifespan.sc,
consumerfrac.sc, endothermfrac.sc, nspp.sc, thermal_bias.sc, npp.sc, human.sc)]
randef <- list(STUDY_ID = ~ temptrend_abs.sc, rarefyID = ~1)
varef <- varPower(-0.5, ~nyrBT)
# full models
if(file.exists('temp/modTfullJbeta.rds')){
modTfullJbeta <- readRDS('temp/modTfullJbeta.rds')
} else {
modTfullJbeta <- lme(Jbetatrend ~ temptrend_abs.sc*REALM +
temptrend_abs.sc*tempave.sc +
temptrend_abs.sc*tempave_metab.sc +
temptrend_abs.sc*seas.sc +
temptrend_abs.sc*microclim.sc +
temptrend_abs.sc*mass.sc +
temptrend_abs.sc*speed.sc +
temptrend_abs.sc*lifespan.sc +
temptrend_abs.sc*consumerfrac.sc +
temptrend_abs.sc*endothermfrac.sc +
temptrend_abs.sc*nspp.sc +
temptrend.sc*thermal_bias.sc +
temptrend_abs.sc*npp.sc +
temptrend_abs.sc*human.sc*REALM,
random = randef, weights = varef, data = trends[i2,], method = 'REML')
saveRDS(modTfullJbeta, file = 'temp/modTfullJbeta.rds')
}
if(file.exists('temp/modTfullHorn.rds')){
modTfullHorn <- readRDS('temp/modTfullHorn.rds')
} else {
modTfullHorn <- lme(Horntrend ~ temptrend_abs.sc*REALM +
temptrend_abs.sc*tempave.sc +
temptrend_abs.sc*tempave_metab.sc +
temptrend_abs.sc*seas.sc +
temptrend_abs.sc*microclim.sc +
temptrend_abs.sc*mass.sc +
temptrend_abs.sc*speed.sc +
temptrend_abs.sc*lifespan.sc +
temptrend_abs.sc*consumerfrac.sc +
temptrend_abs.sc*endothermfrac.sc +
temptrend_abs.sc*nspp.sc +
temptrend.sc*thermal_bias.sc +
temptrend_abs.sc*npp.sc +
temptrend_abs.sc*human.sc*REALM,
random = randef, weights = varef, data = trends[i3,], method = 'REML')
saveRDS(modTfullHorn, file = 'temp/modTfullHorn.rds')
}
summary(modTfullJbeta)
Linear mixed-effects model fit by REML
Data: trends[i2, ]
Random effects:
Formula: ~temptrend_abs.sc | STUDY_ID
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.06301414 (Intr)
temptrend_abs.sc 0.03804902 0.262
Formula: ~1 | rarefyID %in% STUDY_ID
(Intercept) Residual
StdDev: 0.0009169992 0.3174155
Variance function:
Structure: Power of variance covariate
Formula: ~nyrBT
Parameter estimates:
power
-1.350808
Fixed effects: Jbetatrend ~ temptrend_abs.sc * REALM + temptrend_abs.sc * tempave.sc + temptrend_abs.sc * tempave_metab.sc + temptrend_abs.sc * seas.sc + temptrend_abs.sc * microclim.sc + temptrend_abs.sc * mass.sc + temptrend_abs.sc * speed.sc + temptrend_abs.sc * lifespan.sc + temptrend_abs.sc * consumerfrac.sc + temptrend_abs.sc * endothermfrac.sc + temptrend_abs.sc * nspp.sc + temptrend.sc * thermal_bias.sc + temptrend_abs.sc * npp.sc + temptrend_abs.sc * human.sc * REALM
Correlation:
(Intr) tmpt_. REALMMr REALMTr tmpv.s tmpv_. ses.sc mcrcl. mss.sc
temptrend_abs.sc 0.276
REALMMarine -0.918 -0.251
REALMTerrestrial -0.885 -0.259 0.812
tempave.sc -0.005 -0.008 -0.007 0.000
tempave_metab.sc 0.037 0.028 -0.027 -0.038 -0.316
seas.sc -0.025 -0.021 0.036 -0.014 0.206 0.011
microclim.sc -0.008 -0.017 0.010 -0.011 -0.026 0.142 0.149
mass.sc 0.011 -0.008 -0.001 0.001 -0.046 -0.516 -0.055 -0.033
speed.sc 0.014 0.003 -0.023 -0.022 -0.011 0.122 -0.060 -0.015 -0.087
lifespan.sc 0.032 0.042 -0.020 -0.035 -0.020 0.720 0.075 0.034 -0.806
consumerfrac.sc -0.055 -0.055 0.069 0.323 -0.027 -0.080 0.008 0.009 0.050
endothermfrac.sc 0.155 0.082 -0.082 -0.316 0.099 -0.190 0.008 -0.036 0.022
nspp.sc -0.009 0.017 -0.009 -0.019 -0.016 0.027 -0.029 -0.022 -0.162
temptrend.sc -0.006 -0.012 0.006 0.007 -0.099 0.004 -0.031 0.014 0.018
thermal_bias.sc 0.017 -0.010 -0.028 -0.007 0.655 -0.003 -0.093 -0.093 0.030
npp.sc -0.003 -0.027 -0.004 -0.001 -0.216 0.173 -0.237 -0.290 0.009
human.sc -0.115 0.048 0.108 0.101 0.005 -0.006 -0.005 0.003 0.011
temptrend_abs.sc:REALMMarine -0.262 -0.954 0.269 0.244 0.007 -0.021 0.032 0.018 0.011
temptrend_abs.sc:REALMTerrestrial -0.244 -0.859 0.220 0.312 0.011 -0.041 -0.037 -0.019 0.012
temptrend_abs.sc:tempave.sc -0.002 -0.008 -0.001 0.009 0.549 -0.328 0.054 -0.066 -0.084
temptrend_abs.sc:tempave_metab.sc 0.013 0.042 -0.008 -0.019 -0.204 0.724 0.065 0.127 -0.325
temptrend_abs.sc:seas.sc -0.009 -0.037 0.015 -0.016 0.154 0.065 0.715 0.092 -0.049
temptrend_abs.sc:microclim.sc -0.008 -0.025 0.006 -0.009 -0.028 0.135 0.097 0.810 -0.039
temptrend_abs.sc:mass.sc -0.004 -0.008 0.006 0.007 -0.087 -0.361 -0.056 -0.032 0.661
temptrend_abs.sc:speed.sc 0.001 0.008 -0.003 -0.009 -0.025 0.017 -0.014 0.003 0.013
temptrend_abs.sc:lifespan.sc 0.020 0.059 -0.010 -0.018 0.023 0.513 0.048 0.018 -0.553
temptrend_abs.sc:consumerfrac.sc -0.048 -0.084 0.047 0.124 0.032 -0.109 -0.011 -0.001 0.033
spd.sc lfspn. cnsmr. endth. nspp.s tmptr. thrm_. npp.sc hmn.sc
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc 0.194
consumerfrac.sc -0.123 -0.109
endothermfrac.sc 0.065 -0.027 -0.410
nspp.sc -0.062 0.130 0.003 0.040
temptrend.sc -0.004 -0.013 0.006 -0.005 -0.020
thermal_bias.sc 0.028 -0.055 -0.011 -0.015 -0.075 0.001
npp.sc -0.017 0.044 -0.051 -0.062 -0.120 0.024 -0.077
human.sc 0.005 0.005 0.002 0.002 -0.009 0.006 0.008 -0.020
temptrend_abs.sc:REALMMarine -0.008 -0.030 0.058 -0.056 -0.023 0.010 0.009 0.015 -0.045
temptrend_abs.sc:REALMTerrestrial -0.009 -0.035 0.145 -0.166 -0.035 0.008 0.003 0.019 -0.037
temptrend_abs.sc:tempave.sc -0.041 0.012 0.016 0.132 0.020 -0.097 0.036 -0.083 -0.008
temptrend_abs.sc:tempave_metab.sc 0.042 0.477 -0.060 -0.174 0.024 0.012 0.070 0.086 0.001
temptrend_abs.sc:seas.sc -0.024 0.038 0.021 -0.017 -0.015 -0.004 0.087 -0.157 -0.015
temptrend_abs.sc:microclim.sc 0.007 0.017 0.003 -0.047 0.012 0.025 0.002 -0.263 -0.005
temptrend_abs.sc:mass.sc -0.003 -0.569 0.028 0.003 -0.113 0.022 -0.009 0.017 0.014
temptrend_abs.sc:speed.sc 0.363 0.014 -0.053 0.039 -0.047 -0.008 0.004 0.003 0.009
temptrend_abs.sc:lifespan.sc 0.038 0.703 -0.063 -0.040 0.114 -0.010 0.002 0.024 0.004
temptrend_abs.sc:consumerfrac.sc -0.065 -0.100 0.324 -0.130 0.002 0.002 0.001 -0.026 -0.006
tm_.:REALMM tm_.:REALMT tmptrnd_bs.sc:t. t_.:_.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial 0.811
temptrend_abs.sc:tempave.sc 0.006 0.020
temptrend_abs.sc:tempave_metab.sc -0.032 -0.054 -0.498
temptrend_abs.sc:seas.sc 0.050 -0.043 0.169 0.056
temptrend_abs.sc:microclim.sc 0.027 -0.020 -0.043 0.103
temptrend_abs.sc:mass.sc 0.014 0.015 -0.113 -0.475
temptrend_abs.sc:speed.sc -0.022 -0.023 -0.041 0.056
temptrend_abs.sc:lifespan.sc -0.043 -0.053 0.022 0.673
temptrend_abs.sc:consumerfrac.sc 0.096 0.300 0.035 -0.137
tmptrnd_bs.sc:ss. tmptrnd_bs.sc:mc. tmptrnd_bs.sc:ms.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc 0.119
temptrend_abs.sc:mass.sc -0.050 -0.022
temptrend_abs.sc:speed.sc -0.021 -0.006 -0.023
temptrend_abs.sc:lifespan.sc 0.040 0.013 -0.819
temptrend_abs.sc:consumerfrac.sc -0.024 0.002 0.048
tmptrnd_bs.sc:sp. tmptrnd_bs.sc:l. tmptrnd_bs.sc:c.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc 0.088
temptrend_abs.sc:consumerfrac.sc -0.156 -0.148
tmptrnd_bs.sc:nd. tmptrnd_bs.sc:ns. tm.:_. tmptrnd_bs.sc:np.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc
temptrend_abs.sc:consumerfrac.sc
tmptrnd_bs.sc:h. REALMM: REALMT: t_.:REALMM:
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc
temptrend_abs.sc:consumerfrac.sc
[ reached getOption("max.print") -- omitted 9 rows ]
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.6064598 -0.2375081 0.1357510 0.6185365 6.5687599
Number of Observations: 43585
Number of Groups:
STUDY_ID rarefyID %in% STUDY_ID
250 43585
summary(modTfullHorn)
Linear mixed-effects model fit by REML
Data: trends[i3, ]
Random effects:
Formula: ~temptrend_abs.sc | STUDY_ID
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 0.06191799 (Intr)
temptrend_abs.sc 0.03262754 0.297
Formula: ~1 | rarefyID %in% STUDY_ID
(Intercept) Residual
StdDev: 0.004381157 0.2964481
Variance function:
Structure: Power of variance covariate
Formula: ~nyrBT
Parameter estimates:
power
-1.222775
Fixed effects: Horntrend ~ temptrend_abs.sc * REALM + temptrend_abs.sc * tempave.sc + temptrend_abs.sc * tempave_metab.sc + temptrend_abs.sc * seas.sc + temptrend_abs.sc * microclim.sc + temptrend_abs.sc * mass.sc + temptrend_abs.sc * speed.sc + temptrend_abs.sc * lifespan.sc + temptrend_abs.sc * consumerfrac.sc + temptrend_abs.sc * endothermfrac.sc + temptrend_abs.sc * nspp.sc + temptrend.sc * thermal_bias.sc + temptrend_abs.sc * npp.sc + temptrend_abs.sc * human.sc * REALM
Correlation:
(Intr) tmpt_. REALMMr REALMTr tmpv.s tmpv_. ses.sc mcrcl. mss.sc
temptrend_abs.sc 0.285
REALMMarine -0.906 -0.256
REALMTerrestrial -0.874 -0.264 0.785
tempave.sc -0.003 -0.008 -0.008 0.006
tempave_metab.sc 0.047 0.029 -0.033 -0.053 -0.374
seas.sc -0.039 -0.020 0.051 -0.010 0.162 -0.001
microclim.sc -0.010 -0.016 0.012 -0.014 -0.058 0.126 0.119
mass.sc 0.011 -0.011 0.000 0.005 -0.043 -0.498 -0.048 -0.019
speed.sc 0.016 0.001 -0.030 -0.027 -0.019 0.120 -0.070 -0.022 -0.093
lifespan.sc 0.041 0.046 -0.024 -0.044 -0.025 0.698 0.072 0.024 -0.808
consumerfrac.sc -0.060 -0.056 0.056 0.336 0.000 -0.086 -0.006 0.003 0.051
endothermfrac.sc 0.166 0.078 -0.084 -0.340 0.140 -0.213 0.001 -0.038 0.009
nspp.sc -0.012 0.021 -0.008 -0.018 -0.008 0.027 -0.029 -0.021 -0.173
temptrend.sc -0.007 -0.014 0.006 0.009 -0.082 0.003 -0.038 0.018 0.019
thermal_bias.sc 0.026 -0.012 -0.038 -0.010 0.662 -0.032 -0.130 -0.131 0.032
npp.sc 0.000 -0.030 -0.005 -0.003 -0.225 0.183 -0.236 -0.305 0.000
human.sc -0.130 0.059 0.121 0.114 0.008 -0.005 -0.004 0.008 0.012
temptrend_abs.sc:REALMMarine -0.270 -0.949 0.279 0.248 0.005 -0.019 0.028 0.017 0.012
temptrend_abs.sc:REALMTerrestrial -0.249 -0.849 0.221 0.315 0.016 -0.055 -0.041 -0.024 0.017
temptrend_abs.sc:tempave.sc 0.004 -0.006 -0.007 0.012 0.417 -0.273 -0.035 -0.092 -0.081
temptrend_abs.sc:tempave_metab.sc 0.012 0.050 -0.005 -0.025 -0.152 0.623 0.082 0.127 -0.281
temptrend_abs.sc:seas.sc -0.008 -0.040 0.009 -0.023 0.060 0.093 0.561 0.067 -0.039
temptrend_abs.sc:microclim.sc -0.008 -0.029 0.005 -0.013 -0.054 0.128 0.071 0.700 -0.025
temptrend_abs.sc:mass.sc -0.007 -0.013 0.007 0.010 -0.085 -0.305 -0.045 -0.024 0.577
temptrend_abs.sc:speed.sc -0.001 0.011 0.001 -0.007 -0.023 0.007 0.000 0.010 0.025
temptrend_abs.sc:lifespan.sc 0.022 0.073 -0.009 -0.020 0.028 0.434 0.040 0.008 -0.487
temptrend_abs.sc:consumerfrac.sc -0.047 -0.092 0.052 0.102 0.031 -0.112 -0.015 -0.006 0.035
spd.sc lfspn. cnsmr. endth. nspp.s tmptr. thrm_. npp.sc hmn.sc
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc 0.190
consumerfrac.sc -0.132 -0.112
endothermfrac.sc 0.070 -0.020 -0.401
nspp.sc -0.059 0.137 0.008 0.048
temptrend.sc -0.002 -0.015 0.007 -0.007 -0.016
thermal_bias.sc 0.028 -0.061 0.007 0.000 -0.066 0.000
npp.sc -0.003 0.044 -0.030 -0.059 -0.120 0.015 -0.085
human.sc 0.004 0.007 0.005 0.001 -0.011 0.006 0.008 -0.024
temptrend_abs.sc:REALMMarine -0.004 -0.031 0.059 -0.056 -0.023 0.012 0.010 0.020 -0.056
temptrend_abs.sc:REALMTerrestrial -0.007 -0.039 0.129 -0.153 -0.036 0.011 0.002 0.015 -0.044
temptrend_abs.sc:tempave.sc -0.034 0.012 0.024 0.125 0.018 -0.079 0.059 -0.048 -0.009
temptrend_abs.sc:tempave_metab.sc 0.025 0.410 -0.064 -0.163 0.031 0.010 0.055 0.062 0.006
temptrend_abs.sc:seas.sc -0.010 0.029 0.003 -0.043 -0.018 -0.003 0.070 -0.112 -0.021
temptrend_abs.sc:microclim.sc 0.015 0.005 -0.005 -0.052 0.013 0.031 0.001 -0.236 -0.007
temptrend_abs.sc:mass.sc 0.013 -0.499 0.024 -0.002 -0.095 0.022 -0.017 0.019 0.013
temptrend_abs.sc:speed.sc 0.202 -0.007 -0.039 0.029 -0.038 -0.012 -0.001 0.002 0.011
temptrend_abs.sc:lifespan.sc 0.008 0.614 -0.064 -0.032 0.104 -0.010 0.008 0.017 0.007
temptrend_abs.sc:consumerfrac.sc -0.039 -0.105 0.300 -0.094 -0.006 0.005 0.002 -0.011 -0.010
tm_.:REALMM tm_.:REALMT tmptrnd_bs.sc:t. t_.:_.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial 0.795
temptrend_abs.sc:tempave.sc 0.004 0.029
temptrend_abs.sc:tempave_metab.sc -0.036 -0.076 -0.507
temptrend_abs.sc:seas.sc 0.053 -0.061 0.119 0.080
temptrend_abs.sc:microclim.sc 0.033 -0.027 -0.065 0.103
temptrend_abs.sc:mass.sc 0.018 0.028 -0.118 -0.474
temptrend_abs.sc:speed.sc -0.028 -0.029 -0.043 0.060
temptrend_abs.sc:lifespan.sc -0.053 -0.067 0.027 0.663
temptrend_abs.sc:consumerfrac.sc 0.100 0.284 0.041 -0.158
tmptrnd_bs.sc:ss. tmptrnd_bs.sc:mc. tmptrnd_bs.sc:ms.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc 0.116
temptrend_abs.sc:mass.sc -0.044 -0.016
temptrend_abs.sc:speed.sc -0.017 -0.006 -0.027
temptrend_abs.sc:lifespan.sc 0.038 0.007 -0.827
temptrend_abs.sc:consumerfrac.sc -0.032 -0.003 0.058
tmptrnd_bs.sc:sp. tmptrnd_bs.sc:l. tmptrnd_bs.sc:c.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc 0.094
temptrend_abs.sc:consumerfrac.sc -0.164 -0.173
tmptrnd_bs.sc:nd. tmptrnd_bs.sc:ns. tm.:_. tmptrnd_bs.sc:np.
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc
temptrend_abs.sc:consumerfrac.sc
tmptrnd_bs.sc:h. REALMM: REALMT: t_.:REALMM:
temptrend_abs.sc
REALMMarine
REALMTerrestrial
tempave.sc
tempave_metab.sc
seas.sc
microclim.sc
mass.sc
speed.sc
lifespan.sc
consumerfrac.sc
endothermfrac.sc
nspp.sc
temptrend.sc
thermal_bias.sc
npp.sc
human.sc
temptrend_abs.sc:REALMMarine
temptrend_abs.sc:REALMTerrestrial
temptrend_abs.sc:tempave.sc
temptrend_abs.sc:tempave_metab.sc
temptrend_abs.sc:seas.sc
temptrend_abs.sc:microclim.sc
temptrend_abs.sc:mass.sc
temptrend_abs.sc:speed.sc
temptrend_abs.sc:lifespan.sc
temptrend_abs.sc:consumerfrac.sc
[ reached getOption("max.print") -- omitted 9 rows ]
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-6.84526227 -0.37611694 0.04106594 0.54207414 6.16845392
Number of Observations: 42586
Number of Groups:
STUDY_ID rarefyID %in% STUDY_ID
218 42586
coefs2 <- summary(modTfullJbeta)$tTable
coefs3 <- summary(modTfullHorn)$tTable
varstoplot <- unique(c(rownames(coefs2), rownames(coefs3)))
rows1 <- which(!grepl('Intercept', varstoplot) | grepl(':', varstoplot)) # vars to plot in first graph
rows1_2 <- which(rownames(coefs2) %in% varstoplot[rows1]) # rows in coefs2
rows1_3 <- which(rownames(coefs3) %in% varstoplot[rows1]) # rows in coefs3
xlims1 <- range(c(coefs2[rows1_2,1] - coefs2[rows1_2,2],
coefs2[rows1_2,1] + coefs2[rows1_2,2],
coefs3[rows1_3,1] - coefs3[rows1_3,2],
coefs3[rows1_3,1] + coefs3[rows1_3,2]))
cols <- c('black', 'grey') # for Jbeta and Horn models, respectively
offs1 <- 0.1 # offset vertically for the two models
offs2 <- 0.01 # offset vertically for the two models (plot 2)
par(las = 1, mai = c(0.5, 3, 0.1, 0.1))
plot(0,0, col = 'white', xlim=xlims1, ylim = c(1,length(rows1)), yaxt='n', xlab = '', ylab ='')
axis(2, at = length(rows1):1, labels = varstoplot[rows1], cex.axis = 0.7)
abline(v = 0, col = 'grey')
for(i in 1:length(rows1)){
if(varstoplot[rows1[i]] %in% rownames(coefs2)){
x = coefs2[rownames(coefs2) == varstoplot[rows1[i]], 1]
se = coefs2[rownames(coefs2) == varstoplot[rows1[i]], 2]
points(x, length(rows1) + 1 - i + offs1, pch = 16, col = cols[1])
lines(x = c(x-se, x+se), y = c(length(rows1) + 1 - i + offs1, length(rows1) + 1 - i + offs1), col = cols[1])
}
if(varstoplot[rows1[i]] %in% rownames(coefs3)){
x = coefs3[rownames(coefs3) == varstoplot[rows1[i]], 1]
se = coefs3[rownames(coefs3) == varstoplot[rows1[i]], 2]
points(x, length(rows1) + 1 - i - offs1, pch = 16, col = cols[2])
lines(x = c(x-se, x+se), y = c(length(rows1) + 1 - i - offs1, length(rows1) + 1 - i - offs1), col = cols[2])
}
}
NA
NA
Black is for Jaccard total turnover (pres/abs), grey is for Morisita-Horn turnover (considers abundance)
# fig.width = 3, fig.height = 5, out.height=2.5, out.width=3, fig.retina =3 for macbook screen
# double that for external monitor
coefs <- as.data.table(summary(modTfull1)$tTable)
coefs2 <- as.data.table(summary(modTfullJbeta)$tTable)
coefs3 <- as.data.table(summary(modTfullHorn)$tTable)
coefs$mod <- 'Jtu'
coefs2$mod <- 'Jbeta'
coefs3$mod <- 'Horn'
coefs$var <- rownames(summary(modTfull1)$tTable)
coefs2$var <- rownames(summary(modTfullJbeta)$tTable)
coefs3$var <- rownames(summary(modTfullHorn)$tTable)
# extract temperature effects and bind model coefs together
cols <- c('var', 'Value', 'Std.Error', 'mod')
allcoefsfull <- rbind(coefs[grep('temptrend|REALM', var), ..cols],
coefs2[grep('temptrend|REALM', var), ..cols],
coefs3[grep('temptrend|REALM', var), ..cols])
allcoefsfull$var[allcoefsfull$var == 'temptrend_abs.sc'] <- 'temptrend_abs.sc:REALMFreshwater'
# add average temperature effect (across realms) to realm-specific temperature effects
meantempeffect <- allcoefsfull[var == 'temptrend_abs.sc:REALMFreshwater', mean(Value), by = mod]$V1
allcoefsfull[var == 'temptrend_abs.sc:REALMMarine', Value := Value + meantempeffect]
allcoefsfull[var == 'temptrend_abs.sc:REALMTerrestrial', Value := Value + meantempeffect]
# add base temperature:human effect (freshwater) to realm-specific temperature:human effects
allcoefsfull$Value[allcoefsfull$var == 'temptrend_abs.sc:REALMMarine:human.sc'] <- with(allcoefsfull, Value[var == 'temptrend_abs.sc:REALMMarine:human.sc'] + Value[var == 'temptrend_abs.sc:human.sc'])
allcoefsfull$Value[allcoefsfull$var == 'temptrend_abs.sc:REALMTerrestrial:human.sc'] <- with(allcoefsfull, Value[var == 'temptrend_abs.sc:REALMTerrestrial:human.sc'] + Value[var == 'temptrend_abs.sc:human.sc'])
# remove non-temperature effects
allcoefsfull <- allcoefsfull[grepl(':', allcoefsfull$var) & grepl('temptrend', allcoefsfull$var), ]
# add info for plotting
allcoefsfull$lCI <- allcoefsfull$Value - allcoefsfull$Std.Error # lower confidence interval
allcoefsfull$uCI <- allcoefsfull$Value + allcoefsfull$Std.Error
nvar <- nrow(allcoefsfull)/3
allcoefsfull$y <- 1:nvar + rep(c(0, 0.1, 0.2), c(nvar, nvar, nvar)) # y-values
allcoefsfull$varname <- gsub('temptrend_abs.sc:|temptrend.sc:', '', allcoefsfull$var)
allcoefsfull$varname <- gsub('REALM', '', allcoefsfull$varname)
allcoefsfull$varname <- gsub('.sc', '', allcoefsfull$varname)
allcoefsfull$varname <- gsub('^human$', 'Freshwater:human', allcoefsfull$varname)
xlims1 <- c(0.0, 0.05) # for realms
xlims2 <- c(-0.01, 0.015) # for traits
xlims3 <- c(-0.004, 0.0025) # for environment
xlims4 <- c(-0.016, 0.005) # for community
xlims5 <- c(-0.01, 0.015) # for human
ddg <- 0.5 # dodge for each model
set1 <- c('Terrestrial', 'Marine', 'Freshwater')
set2 <- c('mass', 'speed', 'lifespan', 'consumerfrac', 'endothermfrac', 'tempave_metab')
set3 <- c('seas', 'microclim', 'tempave')
set4 <- c('npp', 'nspp', 'thermal_bias')
set5 <- c('Terrestrial:human', 'Marine:human', 'Freshwater:human')
p1 <- ggplot(subset(allcoefsfull, varname %in% set1),
aes(varname, Value, group = mod, color = mod)) +
geom_hline(yintercept = 0, linetype = 'dashed', color = 'light grey') +
geom_errorbar(aes(ymin = lCI, ymax = uCI), width = 0, position = position_dodge(ddg)) +
geom_point(position = position_dodge(ddg)) +
labs(y = 'Temperature change effect', x = '', tag = 'A') +
scale_color_grey() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.position='none',
axis.text=element_text(size=7),
axis.title=element_text(size=7)) +
coord_flip(ylim = xlims1)
p2 <- ggplot(subset(allcoefsfull, varname %in% set2),
aes(varname, Value, group = mod, color = mod)) +
geom_hline(yintercept = 0, linetype = 'dashed', color = 'light grey') +
geom_errorbar(aes(ymin = lCI, ymax = uCI), width = 0, position = position_dodge(ddg)) +
geom_point(position = position_dodge(ddg)) +
labs(y = 'Interaction with temperature change effect', x = '', tag = 'B') +
scale_color_grey() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.position='none',
axis.text=element_text(size=7),
axis.title=element_text(size=7)) +
coord_flip(ylim = xlims2)
p3 <- ggplot(subset(allcoefsfull, varname %in% set3),
aes(varname, Value, group = mod, color = mod)) +
geom_hline(yintercept = 0, linetype = 'dashed', color = 'light grey') +
geom_errorbar(aes(ymin = lCI, ymax = uCI), width = 0, position = position_dodge(ddg)) +
geom_point(position = position_dodge(ddg)) +
labs(y = 'Interaction with temperature change effect', x = '', tag = 'C') +
scale_color_grey() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.position='none',
axis.text=element_text(size=7),
axis.title=element_text(size=7)) +
coord_flip(ylim = xlims3)
p4 <- ggplot(subset(allcoefsfull, varname %in% set4),
aes(varname, Value, group = mod, color = mod)) +
geom_hline(yintercept = 0, linetype = 'dashed', color = 'light grey') +
geom_errorbar(aes(ymin = lCI, ymax = uCI), width = 0, position = position_dodge(ddg)) +
geom_point(position = position_dodge(ddg)) +
labs(y = 'Interaction with temperature change effect', x = '', tag = 'D') +
scale_color_grey() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.position='none',
axis.text=element_text(size=7),
axis.title=element_text(size=7)) +
coord_flip(ylim = xlims4)
p5 <- ggplot(subset(allcoefsfull, varname %in% set5),
aes(varname, Value, group = mod, color = mod)) +
geom_hline(yintercept = 0, linetype = 'dashed', color = 'light grey') +
geom_errorbar(aes(ymin = lCI, ymax = uCI), width = 0, position = position_dodge(ddg)) +
geom_point(position = position_dodge(ddg)) +
labs(y = 'Interaction with temperature change effect', x = '', tag = 'E') +
scale_color_grey() +
theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_blank(), axis.line = element_line(colour = "black"),
legend.position='none',
axis.text=element_text(size=7),
axis.title=element_text(size=7)) +
coord_flip(ylim = xlims5)
grid.arrange(p1, p2, p3, p4, p5, ncol = 2, layout_matrix = rbind(c(1,2), c(3,4), c(5, NA)))